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1. Introduction 

Event-shape distributions in e+e" annihilation are observables which measure the geo- 
metrical properties of energy-momentum flow in a hadronic final state. They have been 
measured over a broad range in energies at LEP [1-4] as well as at earlier electron-positron 
colliders [5,6]. The cvcnt-shapc distributions allow for a detailed probe of the dynamics of 
QCD and especially for a precise determination of the strong coupling constant a^. Owing 
to their infrared and collinear safety, they can be computed systematically in perturbation 
theory. They usually span the range between the kinematical situation of two coUimated 
back-to-back jets (dijet limit) and a perfectly spherical final state. 

The fixed-order description, which expands the distribution in powers of the strong 
coupling constant to leading order (LO), next-to-lcading order (NLO), ncxt-to-ncxt-to- 
leading order (NNLO) and so on, is reliable and convergent over most of the kinematical 
range of the event-shape. In the dijet limit, which is attained for the thrust variable [7] as 
r ^ 1, the convergence of the fixed-order expansion is spoilt by large logarithmic terms 
log(l — T) at each order in the strong coupling constant, which necessitates a resummed 
description. Resummation of the event-shape distribution accounts for the logarithmically 
enhanced terms to all orders in perturbation theory, and ensures a reliable prediction in 
the dijet region. The resummed cross section is organized in terms of leading logarithms 
(LL), next-to-leading logarithms (NLL) and so on. During LEP times, precision studies 
of a standard set of six event-shapes were based on the combination of fixed-order NLO 
calculations [8-15] with NLL resummation [16-18]. To avoid the double counting of terms, 
both expansions need to be matched onto each other and different matching procedures 
are available [19]. 

In the recent past, substantial progress was made both on the fixed-order and the 
resummed description of event-shapes. Following the development of new methods for 
calculations of QCD jet observables at NNLO [20], the NNLO corrections to e'^e" — >■ 3 
jets and related event-shape observables were computed [21-26]. These calculations are 
based on a numerical integration of the relevant three-parton, four-parton and five-parton 
matrix elements, which are combined into a parton-level event generator. Next-to-leading 
order electroweak corrections were also computed very recently [27, 28] . Determinations 
of the strong coupling constant using the newly available NNLO results [29-31] and the 
matched [32] NLL-I-NNLO [33-36] predictions led to a big improvement in the scale uncer- 
tainty and showed the need to go beyond NLL in resummation. 

The resummation of large logarithmic corrections is based on a factorization of the 
event-shape cross section in the dijet limit into a convolution of three leading regions (soft, 
collinear and hard virtual). In the conventional approach [16,37], the resummation of large 
logarithms is accomplished in Mellin (Laplace) space and the resummed distributions are 
obtained by an inverse transformation. In this approach, the NLL corrections to all stan- 
dard event-shape variables were obtained [16-18], as well as NNLL results on the energy- 
energy correlation function [38]. By applying soft-collinear effective theory (SCET, [39]) to 
event-shape distributions [40,41], the resummation can be performed directly in momen- 
tum space. SCET offers moreover a systematic framework to compute all soft, collinear 
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and hard contributions. In this framework, the N'^LL resummation for thrust [42,43] and 
the heavy jet mass [44] have been performed and apphed for a precise determination of 
as, and the framework for the resummation of the jet broadening distributions has been 
outhned [45,46]. In these calculations, the hard subprocess is inferred from the three- 
loop quark form factor [47,48], the collinear jet function is required to the two- loop order, 
which is known from earlier calculations of SCET resummation in heavy meson decays and 
deep inelastic scattering [49,50], while the soft subprocess to two-loop order could be com- 
puted from the renormalization group invariance of the cross section only up to a constant 
term. Using the fixed-order NLO results [15], this term was determined numerically by 
two indepedent groups [42,51], obtaining mutually inconsistent results. Motivated by this 
discrepancy, we perform an analytical calculation of the fiill two-loop soft subprocess for 
thrust from first principles in this paper. A few days before the completion of this work, 
the two-loop hemisphere soft function was derived within SCET [52], and applied to com- 
pute the soft subprocess for the thrust and the heavy jet mass. Shortly after the release 
of our work, two further SCET-based calculations of the logarithmic contributions to the 
two-loop di-jet soft function [53] and of the soft function for the Drell-Yan process [54] 
appeared. 

The major difference between conventional and SCET-based resummation is the han- 
dling of intermediate scales in the calculation. In the SCET-based resummation, these 
scales are fixed to their natural values directly in momentum space, on the other hand in the 
conventional approach they are sampled along a complex contour in peforming the Laplace 
inversion. Although both approaches yield identical results at NLL, owing to the presence 
of a Landau pole in the QCD coupling constant, power-suppressed differences between them 
could appear in higher order logarithmic contributions. These power-suppressed terms are 
outside the scope of the logarithmic resummation, so both approaches could in principle 
yield different but equally correct results. To address the compatibility between the two 
resummations, we perform the NNLL resummation of the thrust distribution in Laplace 
space in this paper. 

The paper is structured as follows: in Section 2, wc review the description of the 
thrust distribution in fixed-order and resummed perturbation theory and establish the 
factorization of the cross section in the dijet limit. Using the eikonal approximation, 
we compute the two-loop soft subprocess in Section 3. The NNLL resummation of the 
thrust distribution in Laplace space and its inversion to momentum space are derived in 
Section 4, and the matching to the fixed-order NNLO distribution is performed in Section 5. 
The numerical impact of these corrections is studied and discussed in Section 6. Finally, 
Section 7 contains our conclusions and an outlook onto future applications. An appendix 
collecting key formulae and technical details of the calculation is enclosed. 

2. Thrust distribution in perturbation theory 

The thrust variable for a hadronic final state in e+e~ annihilation is defined as [7] 




(2.1) 
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where pi denotes the three-momentum of particle i, with the sum running over ah particles. 
The unit vector n is varied to find the thrust direction fiT which maximises the expression 
in parentheses on the right hand side. In the present paper we will mostly work with the 
quantity r = 1 — T. 

It can be seen that a two-particle final state has fixed T = 1, consequently the thrust 

distribution receives its first non-trivial contribution from three-particle final states, which, 
at order ctg, correspond to three-parton final states. Therefore, both theoretically and 
experimentally, the thrust distribution is closely related to three-jet production. 

2.1 Fixed-order and resummed calculations 

The differential thrust distribution in perturbation theory is known at NNLO [22,25]. At 
a centre-of-mass energy Q and for a renormalization scale n takes the form 

^^(t,Q) = a.(M)^(T) + a^(//)^(r,x^) + a3(^)^(r,x^) + 0(a^) , (2.2) 



with 



and where the explicit dependence on the renormalization scale is given by 
dB , . d5_ , . 2>d^. 



-(^'^m) = + 2^0 log(x^) — (r), 

(T,x^)=^(r) + 2 1og(x2)(^2/3o^(r) + 2/3i^(r)) + (2/3olog(x2))^^(r). (2.4) 



For the QCD ^-function we follow the convention given in Appendix A. In theoretical 
computations it is customary to normalize the distributions to the Born cross section 
(To since, for massless quarks, the normalization cancels all clcctrowcak coupling factors. 
However, experimentally it is easier to normalize the distributions to the total hadronic 
cross section a. The correction for the normalization can be done by expanding the ratio 
(To /cr in powers of ctg, which is nowadays known at four loops [55]. In the massless case the 
ratio is given by 



— = 1 + cxsKi + al [K2 + 2/3o log(x^)i^i] + O(a^), (2.5) 



where 



1 

K2 = - 



\cl + CfCa - 44C3) + ^ (-22 + I6C3) 



(2.6) 



The dependence on the renormalization scale is universal and is the same in (2.4) 
and (2.5). Inserting the expansion oi uq/cj we obtain the following expression: 

^^(^'Q) = "^(/^)^W + "'(/^)^(^'^M) + "^(M)^(T,x^) + 0(a^) , (2.7) 
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where A,B,C are related to A, B,C hy 

AiT)=A{T), (2.8) 
B{T,x^)=B{T,x^)-KiAiT), (2.9) 
C{t, x^) = C{t, x^) - K,B{t, x^) + {K^ - K2) A{t) . (2.10) 

For later convenience we consider also the integrated distribution 

R.(r).i/'i!5^e(r-rO<ir', (2.11) 

which has the following fixed-order expansion: 

Rt{t) = l + ^(r)a,(/x2) + B {t, ti^) aUfi^) + C (r, /x^) ^3(^2^^ (2.12) 

The fixed-order coefficients A, B, C can be obtained by integrating the distribution (2.7) 
and imposing i?r(Tmax,Q) = 1 to all orders, where Tmax is the maximal kinematically 
allowed value. 

In the two-jet region the fixed-order thrust distribution is enhanced by large infrared 
logarithms which spoil the convergence of the pcrturbativc scries. The convergence can be 
restored by resumming the logarithms to all orders in the coupling constant. The resummed 
cross section can in general be written as 



Rt{t) = ( 1 + Ckocl J exp \Lgi{asL) + 52(0,^) + — /3o53(a.^) + • • • , (2.13) 
fc=i ^ 



where L = log(l/T). The function gi encodes all the leading logarithms, the function g2 
resums all next-to-leading logarithms and so on. 

The last equation gives a better prediction of the thrust distribution in the two-jet 
region, but fails to describe the multijet region r — )■ Tmax, where non-singular pieces of the 
fixed-order prediction become important. To achieve a reliable description of the observable 
over a broader kinematical range the two expressions (2.12) and (2.13) can be matched, 
taking care of avoiding double counting of logarithms appearing in both expressions. To 
this purpose we reexpand the functions gi in powers of a^L: 

Lgi{asL) =Gi2asL^ + G2^alL^ + Gs^a^^L^ + ... , 
g2{asL) =GiiasL + G22a^L^ + Gssa^L'^ + ... , 
^PogsiasL) =G2ia^L + Gssa^L^ + G^sajL"^ + ... . (2.14) 

TT 

In a second step we reexpand the exponential function recovering the full logarithmic 
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dependence of the fixed-order series at NNLO 

^log (^) =Ci + GuL + GuL^ , (2.15) 

+ (G23 + G12G11) + ^G?2^^ , (2.16) 

R^g (r) =C3 + (G31 + C1G21 + C2G11) L 

+ ^^32 + C1G22 + 2^1^11 + C'2Gi2 + GiiG21^ 
+ ^^33 + G11G22 + G12G21 + G1G11G12 + Q^ll + G1G23 ) L"^ 



+ ( G34 + G12G22 + nCiG\2 + G11G23 + -^G\iGi2 I L'^ 



-CIG12 + G11G23 + -( 

+ (G12G23 + ^G^^Gii) + ^Gl2L^ . (2.17) 

The difference between the logarithmic part and the fuh fixed-order series at different orders 
is given by 

di(r)=^(r)-<)(r), (2.18) 
d2(r)=H(r)-i?;2(r), (2.19) 
d3(r) =C(r) - i2f2(r) . (2.20) 

The functions di(r) contain the non- logarithmic part of the fixed-order contribution and 
vanish for r — > 0. We collect them into a function D{t, Q) defined as 

D{t, Q) = asdi{T) + a%{T) + a%{T) + C>(a^). (2.21) 

2.2 Kinematics and factorization of thrust 

The definition of thrust (2.1) splits the final-state into two hemispheres Sfi^. and Sfi^. 
separated by the plane Pfj orthogonal to the unit vector Ht- Each final state particle with 
momentum pi is assigned to either Sfi^, or S^j^, depending on whether pi ■ nx is positive or 
negative. As it was shown in [16], no final state momenta can lie in P^. Wc denote with p 
and p the total momenta in the hemispheres Sfi^ and respectively. We can parametrize 
the two total momenta as 

where n'^ = {I^ut) and = (1,— n-r)- A simple kinematical analysis [16] shows that in 
the dijet limit we can recast (2.1) as 

r^l-T=^ + ^ + 0{{1 - Tf), (2.23) 
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where we are neglecting terms of relative order 1 — T which give rise to power-suppressed 
contributions to the cross section. In our analysis we need to separate soft particle contri- 
butions from collinear ones, so we have to modify (2.23) in order to single out the explicit 
contributions arising from each of the two kinematical configurations. Let us then con- 
sider a hard parton with momentum in the Sfi^ hemisphere which produces a hard 
collinear final state parton with momentum k^^ after emitting a soft gluon with momentum 
(Fig. 1). The parent parton is moving along so its hard momentum component p ■ n 
is of order 0{Q), while the remaining components p ■ n and \p±\ are suppressed. The soft 
gluon momentum components tend to zero with the same scaling (q ■ n q ■ n ^ \q±\). It 
is then easy to see that ^ p'^ — Qq • n. Plugging the expression for p'^ arising from the 
previous equation into (2.23), we end up with the expression 



^ k'^ k'^ q ■ n q ■ n ^. 
1-T = \ h h h 0((l - T) 



2^ 



where k (k) is the total collinear momentum in the S^^, {Sa^) hemisphere and q (q) is the 
total soft momentum in the S'^j, (Sfij,) hemisphere. 




Figure 1: First branching kinematics. 



An important simplification in the calculation is achieved if the recoil effects due to 
emissions of wide angle soft particles are neglected. Before neglecting them we have to 
make sure that they do not give rise to any logarithmically enhanced terms. To see it, we 
rewrite the thrust as 

l-^ = E^(l-^o^^*)' (2.25) 
i ^ 

where Wj is the energy of the i—th final state particle and 9i is the angle of its direction to 
the thrust axis. We now consider the wide angle soft contributions to T, but we neglect 
their effect on the determination of the thrust axis ut- We call this fake thrust axis nfake 
and we define 5 as the small angle of nfako to the physical axis ut- In approximating ut 
with nfake the angles 9i are replaced by 9[ = 9i — 5. 

We now consider the angle 6 due to a single large angle soft emission. From simple 
kinematics we obtain 6 ~ oJs/Q, where ajg is the soft particle energy. We observe that 6 
is of the same order as the soft emission contribution to the thrust (i.e.~ r). We want to 
estimate the effect of neglecting the recoil 5 on the thrust itself using expression (2.25). The 
effect of the approximation on collinear emissions {9 ~ 0, cjc ~ Q) is At ~ {uJc/Q)S'^ ~ r^, 
while for large angle soft emissions {9 » 0, a;^ <C Q) we find At ~ {uJs/Q)6 ~ r^. 
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We see that in both cases the effect of the recoil amounts to a contribution to T of 
relative order (at least) 0(1 — T), so it produces power-suppressed terms. We can then 
replace the thrust axis in the dijet region with the direction of the hardest (jet-initiating) 
parton. From now on this approximation is understood. 

Factorization properties of event-shapes have been widely studied in the literature [40, 
56,57]. Referring to Fig. 2 we recast the cross section (2.11) as 

Rt{t)= H(^,as{fi)] [ dk^dPj{-,as{fi))Ji-,as{fi)) 
VA' J J fJ- 

X [ dwS{-,asip))e(Q^r -k^ -k^ - wQ) + Oir). (2.26) 

J ^J' 

We use the integral representation of the 0-function 

G(qV -p-k^-wQ) = — f — e-<3%-fc%-S%-""9, (2.27) 

27r^ Jc V 

and the Laplace transform to recast Eq. (2.26) as 

where we set = uQ^ and A^o = e~''^ . For the sake of simplicity we defined CTAr((5^, a<j) 
as the Laplace-transformed cross section and we omitted the term 0{t). The soft subpro- 




Figure 2: Leading regions in dijet factorization. 

cess S {Nq/NQ/ fi, UsiiJ,)) describes the interaction between the two jets of hard collinear 
particles through soft gluon exchange. It can be therefore defined in a gauge invariant way 
as a correlator of Wilson lines 

S(^^,as{fi)] = ^ I dTse--^^Y.^mkml{me^k)Jcut{TsQ){keik\Wn{mn{m). 

(2.29) 
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where we defined Ts = w/Q. Wn and Wn are Wilson lines 

Wn{y) = Pexp (^ig dsn- A{ns + y)^ , (2.30) 

describing the eikonal interaction of soft gluons with the fast moving quarks along the 
directions and respectively. A{ns + y) in the previous expression denotes the gluon 
field in QCD. The sum runs over the final states \keik) involving k soft particles whose 
phase space is constrained according to the thrust measurement function Jcut{'''Q'^) ■ Both 
soft and soft-coUincar contributions are encoded into the soft subprocess. The one-loop 
expression of the soft subprocess is known since a long time and we compute it with two- 
loop accuracy in the next section. 

The collinear subprocess J {J) describes the decay of the jet-initiating hard quark 
(antiquark) into a jet of collinear particles moving along the (n^) direction. It is 
therefore an inclusive quantity which can be found in many other relevant QCD processes 
such as deep inelastic scattering and heavy quarks decay [37,49,58]. Double-counting of 
soft-collincar and n-coUincar (n-collinear) contributions has to be avoided when computing 
the n-coUincar (n-coUincar) jet subprocess. To this end different regularization schemes can 
be found in the literature and they all provide the same results for the logarithmic structure 
of the jet subprocess. For the purposes of the present paper the explicit expression of the 
collinear subprocess is not required. All we need is its non-logarithmic term at two-loop 
order, which has been computed in [49] using dimensional regularization. In that work, all 
pure virtual corrections to the collinear subprocess are given by scaleless integrals and thus 
they vanish in dimensional regularization. As we will see in detail in Section 3 this property 
holds true also for the soft subprocess and it ensures that the whole virtual contribution is 
encoded into the hard subprocess defined below, which can be identified with the squared 
of the constant part of the quark form factor. 

The short-distance hard function H {Q/fi,as{n)) = IH {Q / n, as{n)) P takes into ac- 
count the hard virtual corrections to the quark-antiquark production subprocess. It is free 
of large logarithms and it can be generally defined such that Eq. (2.26) reproduces the 
fixed-order cross section up to power suppressed terms. Since the hard function contains 
only constants, one can avoid it by performing a log(i2)-matching to fixed-order as it will 
be shown in Section 5. We nevertheless take it into account in order to compute the full 
constant part of the cross section with 0{a^) accuracy. 

3. Soft gluon emission at two- loop level 

Let us turn into the computation of the soft subprocess (2.29). Expanding the generating 
functional to the desired order we obtain the relevant set of cut Feynman diagrams. For 
real emissions, the integration runs over the soft gluons phase space constrained by the 
thrust measurement function in the dijet limit. The leading order contribution trivially 
reduces to S{ts). In what follows we compute S{Ts,Q/n,as{iJ,)) to two-loop level. 
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3.1 One-loop result 

At one- loop level the evaluation of the Wilson loop is straightforward. The contributing 
diagrams are listed in Figure 3. 



Figure 3: NLO contribution to the soft subprocess. 



The virtual soft correction (Fig. 3(a)) to the vertex identically vanishes in dimensional 
regularization since the integration over the loop momentum is given by the scaleless inte- 
gral 



(27r)'^ {q-n + iO+){q ■ n - iO+){q^ + ^0^ 



-, (3.1) 



where the phase space constraint reduces to 5{ts) since no real gluons have been emitted. 
The real contribution involves the emission of a soft gluon off one of the two Wilson lines. 
The emitted gluon may go either into the S^^, hemisphere or into the 5^^, one. The resulting 
phase space measure is then 



d'^q 



(27r)5W(g2)(5(r,Q - q • n)e{q ■ n - q ■ n) + 5{tsQ - q • n)e{q -n-q- n)), (3.2) 



where the two ©-functions forbid the emitted gluon from going backwards, heading towards 
the opposite hemisphere. The integrand is the eikonal factor corresponding to the single 
emission shown in Figure 3(b) and its mirror conjugate diagram. The sum of the two then 
reads 



^{tsQ — Q ■ n)Q{q ■ n — q • n) + 5{tsQ — q ■ n)@{q ■ n — q - n) 
{q ■ n + iQ+){q ■ n + iQ+) 



We first replace 



d" 



-d{q ■ n)d{q ■ n)d'^ "^q^ 



(3.3) 



(3.4) 



and evaluate the integral over Since there is no explicit dependence on q± in the 
integrand function, it simply reduces to the replacement 



d''-\^5^+\q^) 



vr 



l-e 



r(i-6) 



{q-n) ^{q-n) ^Q{q ■ n)@{q ■ n). 



(3.5) 
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It is then straightforward to show that the final result reads 



/X / TT er(l — e) \fi 



asifJ.) 



Q 



-2e 



l+2e 



(3.6) 



where we set n-n = 2 and we replaced the bare coupling = g'^/An with the renormalized 
one in the MS scheme, i.e. 



2e 



(4vr) 



11Ca-2Nf asiii) . > 



with 7b = 0.5772 . . . being the Euler-Mascheroni constant. Using the relation 

l+2e 



2e 





"1" 


-2e 


"logr/ 


+ 








Js. 


+ 





we recast (3.6) as 

S«(r..2,„,(,)).c.^^(_^ + ? 



*(r.)(^-21og=i^ 



+ 0(e'), 



log^ + logrs 



-1 + 



The above result in Laplace space reads 



S^' (^,a.(/^)) = j^^'cir.e-^-^W (r.,^,a.(M) 



Cf 



21og^^Y 
4 ^ NoQj 



(3.7) 



(3.8) 



(3.9) 



(3.10) 



where ATq = e"'^^ . 

Note the absence of single logarithms in the 0{e^) term meaning that the one loop 
soft subproccss receives logarithmic contributions only when the emitted gluon is both soft 
and collincar to one of the eikonal legs. Wide angle soft emissions do not contribute at this 
order, but they become relevant at two-loop level. 



3.2 Two-loop result 

At order ©(a^) up to two real gluons are emitted. Note that diagrams with virtual dressing 
of eikonal lines vanish identically since the Wilson path lies on the light cone {i.e. = 
fp- = 0). Furthermore, we only need to evaluate those diagrams contributing with maximal 
non-abelian (i.e. CpCa) and fermionic (i.e. CpTpN^) terms, because of the non-abelian 
exponentiation theorem [59,60]. According to it we can write 
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where s^'^\NQQ/{Nn),as{iJ-)) involves CfCa and CpTpNp contributions while Cp terms 
arise from exponentiation of the 0{as) result S^^\NoQ/{N^),as{^J-)). It follows that the 
C(a^) result is given by 

where S^^^NoQ/{Nn),as{ij)) was computed in the previous section. The remaining set 
of (non-vanishing) diagrams contributing to s^'^\NoQ/{Nfj,),as{fJ.)) is shown in Figure 4, 
where the mirror conjugate diagrams are omitted. The vacuum polarization blob includes 
fermions, gluons and ghosts as depicted in Figure 5, since the calculation will be carried 
out in the Feynman gauge. 




(b) (c) (d) 




(f) (g) (h) 




(i) (j) (k) (1) 

Figure 4: NNLO contribution to the soft subprocess. Grey blobs stand for the sum of vacuum 
polarization bubbles due to fermions, gluons and ghosts. To complete the set one has to take into 
account mirror conjugate diagrams in addition. Pure abelian diagrams (i.e. proportional to Cp) 
are omitted. 





Figure 5: Vacuum polarization diagrams in the Feynman gauge. 



The double virtual contribution (Fig. 4(i-l)) is made of scaleless integrals which iden- 
tically vanish in dimensional regularization. 

The first non-trivial contribution to consider is the one-loop virtual diagram with an extra 
real gluon, depicted in Fig. 4(e,f). The sum of such diagrams and their mirror conjugate 
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ones is gauge invariant. Using the known expression of the soft current at one loop [61], 
we can write their contribution as 



(2) ( Q r ^\ ,9^1, A^er^r. I {I - e)T^ {I + () 

(f\ ^s, -,as{l^) = -4.^{47rfiyQCFCA ^ — - 



X 



/ 



d'^q (27r)5W(g2) (5(rsQ - q ■ n)Q{q ■ n - q ■ n) + 6{tsQ - q ■ n)&{q -n-q-n)) 
(2^ {{q-n){q-n)Y+^ ' 

(3.13) 



where the phase space constraint is the same as in the one loop case, since only one real soft 
gluon is emitted. Notice that the pole prescription (±iO) has been omitted since the poles 
are never "touched" during the integration because of the B-functions in the numerator. 
Prom now on we write explicitly the factor (iiO) only when relevant. Using the relation 
(3.4) and evaluating the straightforward integral left we obtain 

(2) / Q ,,\_ asilJ^fCFCAe'^^^ T\l-e)T^{l + e) ( i^^' ( 
^e),(f)\^^^ ^^"^sWj - ^2 2 e3 r2(l-2e)r(l + 2e) ' 

(3.14) 



where we expressed the bare coupling according to equation (3.7). 

Let us now consider the double real emission. The phase space constraint is more 
involved since two distinct gluons are emitted. Each of them can indeed go either into the 
Sfi^ or the Sfi^ hemisphere, leading to the following phase space cut 

Q{2'K)5^+\q^){2TT)5^+\e)Jcut{TsQ) = g(27r)5(+)(g2)(27r)5(+)(A:2) 

X {5{tsQ — q ■ n — k • n)@{q ■ ft — q ■ n)Q{k ■ n — k ■ n) (3.15) 

+ 5{tsQ — q ■ n — k ■ n)Q{q ■ n — q ■ n)@{k ■ n — k ■ n) (3.16) 

+ S{tsQ — q ■ n — k ■ n)Q{q ■ n — q ■ n)Q{k ■ n — k ■ n) (3-17) 

+ 5{tsQ — k ■ n — q ■ fi)Q{k ■ n — k ■ n)Q{q ■ n — q - fi)). (3.18) 

Such a constraint gives rise to non-trivial phase space integrals which are solvable using 
analytic techniques. We evaluate the integrals in MATHEMATICA, partly using the package 
HypExp [62]. An independent numerical check is also performed using the computer code 
SecDec [63], a recent implementation of sector decomposition based on the algorithms pre- 
sented in [64,65]. The numerical integration is then carried out using the Monte Carlo 
routines BASES [66] and VEGAS [67] included in [63] with 10^ Monte Carlo events per coef- 
ficient. The two results agree within an uncertainty of 0.001%. 

We organise the calculation considering first the class of diagrams without any internal 
gluons {i.e. gluons which arc not involved in the cut) shown in Fig. 4(a,b) and mirror 
conjugate diagrams, then the class of diagrams with only one internal gluon (Fig. 4(c,d) 
and mirror conjugate diagrams) and finally vacuum polarization diagrams, containing two 
internal gluons (Fig. 4(g,h) and mirror conjugate diagrams). 
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3.2.1 Box-type diagrams 

We first consider the Box- type diagrams depicted in Fig. 4(a,b). The second diagram 
(Fig. 4(b)) has a simple structure and can be evaluated easily. Dropping the abelian part 
of the colour factor we are left with 



d'^q d^k Q{27r)5(+Hq^){27T)6M{e)Jcut{TsQ) , . 

(27r)'^ (27r)<i {q ■ n){q ■ n){k ■ n){k ■ n) ' ^ ' 

where Jcut{TsQ) is the sum of terms in round brackets defined in Eqs. (3.15-3.18). The in- 
tegrand function docs not depend on the transverse component of the integrated momenta, 
so we can use Eqs. (3.4,3.5) getting 

(2) / Q , A CaCf^ {n-nf 4 4e ^'"'^ 
-(4^--'«^(/^)J= 2-^(2^^^ 4r2(l-e) 

X / d{q ■ n)d{q ■ n)d{k ■ n)d{k ■ n 



S{tsQ 


- q 


■ n — k ■ n)&{q ■ n 


-q-n)@{k 


■ n — 


k' 


■ n) 
—+ 




(q- 


ny+^{q ■ ny+^{k ■ 


■ n)^+'^(/c • n] 


,1+. 






6{tsQ 


- q 


■ n — k • n)Q{q ■ n 


— q ■ n)Q{k 


• n — 


k' 


■ n) 
— + 




(q- 


n)i+'=(g • n)i+^(fc • 


■ n)i+'^(fc • n] 


,1+6 






S{tsQ 


- q 


■ n — k ■ n)Q{q • n 


— q • n)Q{k 


• n — 


k- 


■n) 



SijsQ — k ■ n — q ■ n)@{k ■ n — k ■ n)G{q ■ n — q ■ n) 
{q ■ ny+^{q ■ n^+^ik ■ n^+^ik ■ n)i+^ 

(3.20) 

We now analyse each of the integrals in Eq. (3.20). We first integrate out the k-n component 
in the first and fourth integral and the k ■ n component in the second and third integral by 
using the (5-functions. We then make the following four changes of variables 

k • n Qtsu{1 — t) q - n ^ QTgts q • n ^ QTgt, (3-21) 

k ■ n ^ Qtsu{1 — t) q • n ^ Qrgts q • n ^ QTgt, (3.22) 

k ■ n ^ Qtsu{1 — t) q- h ^ QTsts q • n ^ QtsI, (3.23) 

k ■ n ^ Qtsu{1 — t) q • n ^ QTsts q • n ^ QTgt, (3.24) 

in the first, second, third and fourth integrals of Eq. (3.20) respectively, followed by the 
replacements s — > ^ , u ^ Mn each of them. We end up with four identical integrals on a 
three dimensional unit cube. Summing them up we obtain the following expression 



"1 



AQtsY dtduds{l-t)-'~^H-'-^'u-'+'s-'+', (3.25) 

Jo 

that can be easily evaluated in terms of Euler beta functions yielding 

{2) f Q I \ CACp^in-nf 4 4, ^2-2e i_4, 1 

ib) [^^^ -' «^(^) ) = —^ {27ryd-2 9 f r2(i_,) (Q^s) -^Bi-2e, -2e). 

(3.26) 
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After replacing the bare coupling with the renormalized one (Eq. (3.7)) and some algebra, 
we can recast the previous equation as follows 

-(4^--'«3(m)J=-^ Vs) e^V^{l-e) ' ('-'^^ 

where we set n • n = 2. 

We now consider the box diagram in Fig. 4(a). Wc note that the mirror symmetrical 
diagram has the same expression provided wc exchange k -H^ q, so it can be taken into 
account by just including a factor of two. The full result reads 



/ 



dS d'^k Q(27r)<5(+)(g2)(27r)<5(+)(fe2)J,„i(r,Q) 



(3.28) 



(27r)'^ (27r)'^ {{q + A;) • 7l){q ■ n){q ■ n){{k + q) ■ n) 
Furthermore, exploiting the symmetry of the integrand under the transformation 

{A;-n-f->-A;-n, g-nOg - n}, 

we see that the integrals arising from the terms (3.15) and (3.17) equal those arising from 
(3.16) and (3.18) respectively. Using the parametrization shown in Eq. (3.21-3.24), we are 
led to the following expression 

(2) / Q , A CaCf asifif ffiY'fl\ e^T^^ 



^^^ys, ^,^sy^jj- 2 ^2 yQj y^j r2(i-e) 

f f\ , , (1 - t)-^-^H^-^'s-^+'u' f\ , , (1 - t)-^-^H^-^'s-^+'u' \ 

X / dtduds — r h / dtduds-. 7- ttt- 7- rr • 

Vio tu + s{l-t) Jo {t + s{l-t)){l-t{l-u))J 

(3.29) 

The two integrals in Eq. (3.29) can be easily evaluated with the desired accuracy, yielding 

/ dtduds^- ; ^ = — -+4C3 + -^e + , 3.30 

Jq tu + s{l — t) be 9 



dt du ds Y Jr. L ^^ = A + + SCa + + 0{e^). (3.31) 

{t + s{l - t)){l - t{l - u)) 2e3 3 e 5 ^ ^ ^ ^ 

Plugging them back into Eq. (3.29) we find 

2^ + y^ + 10C3 + ^7r^6 + O(e^)). (3.32) 

This completes the evaluation of Box-type diagrams contribution to the two loop soft 
subprocess. 
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3.2.2 Non-abelian diagrams 

The class of non-abelian diagrams (Fig. 4(c,d) and mirror symmetrical diagrams) is much 

more involved due to the presence of the three gluon vertex. The mirror conjugate diagrams 
in which the real gluon is connected to the opposite leg to the right of the cut are related 
to those depicted in Fig. 4(c,d) by the transformation 

{k-n-(-^q-n,k-n-^q - n}. 

They can be taken into account by including a factor of two. 

Thus we can write the whole non-abelian contribution as follows 



(2) 

s) 



(cW) yrs,-,as{^^)j=C,CA^i g {n ■ n) J (2vr)^-i (2vr)^-i « ((, + fc)2 _ ,o) 

2k-n + q-n q-n — k-n 

xJcut{TsQ){^j^ • n)(g • n){{q + k) ■ n) ^ {{k + q) ■ n){k ■ n){{k + q) ■ n)'^ 

{k ■ n ■H' q ■ n,k ■ n ■H' q ■ n}). (3.33) 

Looking at the expression of the phase space constraint Jcut we see that the transformation 

{k ■ n q ■ n, k ■ n q ■ n} 

has only the effect of exchanging the terms (3.15) and (3.17) with (3.16) and (3.18) re- 
spectively, so the contribution due to the last term in round brackets amounts to multiply 
once again by two. Unlike the box-type case, the integrand function depends explicitly on 
the transverse components through the gluon propagator. We use (3.4) and perform the 
integral over the transverse components as shown in Eq. (C.6). Writing the result in terms 
of the parametrization (3.21-3.24) we end up with the following expression 

'icUd)[Js, ^,asWJ - 2 5 (27r)2d-2 l^gj |^^J r(l - e) r(l - 2e) 

X l\tduds( (1 - t)~^-H-^-^^s-^+^u-^+^{s + u) ~ ~ + - ^) 
Jo V s{l-t) + tu 

2F1 (ll-e,l- 2e, 4^7^^) 



j(l -t){2- t{2 - u)) + t{l - t{l - 2u)) 2-^1 (l' 5 - 1 - 2e, 4 



{s{l-t) + t){l-t{l-u)) {\ + ^f^f 



(3.34) 



where the first integral is made of the contributions due to the terms (3.15) and (3.16), 
while the ones due to (3.17) and (3.18) are encoded in the second integral. Since s < 1 and 
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■u < 1 we can use Eq. (C.IO) to perform the following replacements 



2Fi(l i-6,l-2e,4 ^ 



F, 1 



2-f 1 



= {^^s + V^f (ois - u)^ 2Fi{l, 1 + e, 1 - e, ^) + 

1 s \ 

+ e(n - s)- 2F1 (1, 1 + e, 1 - e, -) , (3.35) 
u u J 

(^1, ^ - e, 1 - 2e,4-^^-^^^=pj = (1 + V^f2Fi{l, l + e,l-e,su), 

(3.36) 



in the first and second integral of (3.34) respectively. 

Identity (3.35) splits the first integral in (3.34) into two simpler integrals. Wc substitute 
It — )• in the first of such integrals and s — )• in the second one in order to remap the 
integration range to the unit cube and we sum them up. After renaming the variables we 
finally obtain 

%Ud) \ ^s, -,as{fi) ] = -^{n • n)- 



m' 7 8 ' ' 7r2 r(l-e)r(l-2e) VQ 



Jo V ' ^ ' {l-t{l-z)){z + t{l-z)) 

X2Fi(l, 1 + e, 1 - e, z) + (1 - t)-^-^H-^-^'u-^+'z-^+'{l + uz) 
z{l - t){2 - t{2 - u)) + t{l - t{l - 2u)) 



{z{l-t) + t){l-t{l-u)) 



-2Fi(l,l + e,l-e,2;tx)y 



(3.37) 



where we used Eq. (3.7) for the coupling. The two integrals in round brackets of Eq. (3.37) 
are given by 

7o V «^ « " ^ ^^ + ^> ^i_t{l-z)){z + t{l-z)) 

A 1 r(l-e) / I 2 4 + 7r2 

- i (6 + + 6C3) - ^ {Ui + 2<1t' + tt" + I8OC3) e + 0(e^)) , (3.38) 

r(l -e) f 1 1 f,,^ 2_2 



X 2F, (1, 1 + M - e, zu) J = ; ^^^^^^^^ ( - Ti^ - 7 - r 

- ^ (207r2 + 97r^ - 42OC3) + C(e)^ . (3.39) 
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Plugging them into (3.37) and setting n • n = 2 we find 



+ I 8 - TT^ + 1 + 



IOtt^ 1997r^ 124C3\ 
3.2.3 Vacuum polarization diagrams 

The last class of diagrams we have to take into account to complete the computation of 
the double real radiation contribution involves diagrams of the type shown in Fig. 4(g,h). 
Summing up diagrams (g) and (h) and their mirror symmetrical ones we end up with the 
following expression 

{{q • n){k ■n)-{k- n){q ■ n)f 2(1 - e)C7^CF - ACpTpUF 



+ i 



{{k + q) ■ nf{{k + q) ■ n)^ {{k + qf + iQ){{k + qY - iO) 
1 CFTFnF-CACF\ 



{{k + q) ■ n){{k + q) ■ n) {k + q^+iO J' 

(3.41) 



where the phase space constraint J^cutiQTs) is the usual measurement function defined in 



(3.15-3.18). We name the two integrals appearing in (3.41) //qW^^] and -^S'S respectively 



and we evaluate them below. We first use (3.4) getting 

j(2),{a) ^1 f d{q ■ n)d{q ■ n) d{q ■ n)d{q ■ n) {{q ■ n){k ■ n) - {k ■ n){q ■ n)f 



'i9),{h) 4 7 (27r)<i-i (27r)<^-i Q'^'^tiQ^s) + ^) . ^)2((^ + . 



{{k + qf + iO){ik + qf - iO) ' 

r(2),(b) _ 1 f d{q- n)d{q ■ n) d{q ■ n)d{q -n) ^^ 
\9),ih) -4 J (27r)^-i (27r)^-i ^-^cutWr, 



I 



{{k + q) ■ n){{k + q) ■ n) 
'<'-\±d''-^q±6M{q'^)S(+){k'') 



(3.42) 



(3.43) 



(A; + q)^ + iO 

From the symmetry of the previous expressions under the transformation 

{k-n-(r^k-n,q-n-(r^q - n}, 

we see that the terms arising from the constraints (3.16) and (3.18) are identical to those due 
to (3.15) and (3.17) respectively. The two internal integrals over the transverse components 
of the soft gluon momenta are evaluated in (C.9) and (C.6). Following the same technique 
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used for the non-abelian integrals and using the following two identities (obtained from 
Eq. (C.12)) 



2^1(2,1-6,1-26,4 



2 --(Vi+v^)0=^^+^^'(®^^-")^^''^('''+^''-^'3+ 

+ @{u - s)ir2Fi f2, 1 + e, 1 - e, - 



(3.44) 

Fi (^2, 1 - e, 1 - 2e, 4^^-^^^=^^ = (1 + ^/^)SFl (2, 1 + e, 1 - e, su) , (3.45) 



in addition to (3.35,3.36), after renaming some of the variables we recast as 



(2),(„) _ 1 VT^- vri- 4- r^(i-e) / 1 ^^1^ / 1 N 1+4^ 



2 r(l - e) r(i - e) (271)2^^-2 r(l -2e)\Qj \ts 

X r dt du dz (2(1 - ^)l-2e^l-2e^-l+26^e (1 - ^j^^ (2, 2 + £, 1 - e, z) + 

(1 - zuf 



+ (1 - t)^-^¥-^W . .^11,2^^1^ ,.^^2 2-?^l (2, 2 + 6, 1 - 6, .^) ) , (3.46) 



and as 



ji2m 1 ^^-^ ^^-^ 4- r2(i-e) ^^^Vl^'^'^ 



2 r(l - e) r(i - e) (27r)2^-2 r(l 2e) \Q J \t, 

rl ( (1- /)-2^y--2e„-l+2e e 

X / dtdudz[2 ^^ n % Fi(l,l + 6,l-6,z) + 

7o V t + 2:(l-t) 

f 1 - t)~'^H~'^^u^z^ \ 

+ (t+;(i-t))(i-«(i-u)) '^'(^''+'^-'-'^' 
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For the four integrals appearing above we obtain 

/ dtdudz^ 2i^i(l,l + e,l-e,^;) - ^ 

Jo + -L. -yl \ — t\ 



t + z{l-t) ^ ' ' ' ^ 2er(l + e)r(-2e) 

^(^^l^(^^i)^h^4^^Uo(.')l (3.48) 



V4e2 2e 

\ r , (1 - t)-^H-^'u'z' ^ , . r(l - e) 

dtdudz -7- :-r2Fi{l,l + e,l-e,zu) - 



{t + z{l-t)){\-t{l-u)y ' ' ' ' ' r(l + e)r(-2e) 

I dt du dz{l - ty-''H^-^'u-^+^'z' . . ^ 2F1 (2, 2 + e, 1 - e, z) = 

1 r(l-e) /I 5 16 7r2 



Jo 



2er(2 + e)r(-l-2e) \l2e^ 36e 27 18 

59 57r2 5C3 

^ h — 

,81 54 6 

dt dudz{l - tf-H^-^'z'u'-- TT2 2^1 (2, 2 + e, 1 - e, z^^) 



+ + + + (3.50) 



{z{l-t) + tY{l-t{l-u)f 
r(l-e) /I /I 7r2\ 7 57r2 7C: 



r(2 + .)r(-i - 2.) V . V6 T8 J - 9 - ^ + f + "f^' J • 

We finally plug these expressions back into (3.41) and we obtain the following result for 
the vacuum polarization diagrams 

4e / 1 \ l+4e 



^ ^ 1 5 /28 7r2\ /20 377r2 62C3\ > 

5 , 31 , /47 57^2^^ , f211 , 1557r2 155C3 \ 2^ 



(3.52) 

where we used (3.7) to replace the bare coupling with the renormalized one in the MS 
scheme. This completes the computation of the relevant contributions to the two-loop soft 
subprocess. As a further eheek of our calculation, we observe that summing up all the 
integrand functions contributing to the two-loop soft subprocess we reproduce the known 
double-soft current derived in [68]. 

3.3 Renormalization of the two-loop soft subprocess 

Before the subtraction of the overall divergences we need to handle the subdivergences. 
The coupling renormalization (3.7) leads to a counter-term 

(2) / Q , A OLs{ii) jio (1) ( Q I ^ "s(^) 

«c.t. —,as{lJ.) = 5"^ ' Ts, —,as{jl) = -Cf 2~ 

V /V, J TT e \ J TT^ 



7r2 e2r(l-e) 

, , , , ^Ca-^Tfuf], (3.53) 
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where S^^'^ (ts, Q / fJ,, as{iJ-)) is the one-loop contribution (3.6). 

After performing the Laplace transform of Eqs. (3.14), (3.27), (3.32), (3.40), (3.52), (3. 53) 
by means of the relation 



e-'-r(-..,(||)", (3.54) 

and summing them up we obtain the following expression for the non-abelian part of the 
unrenormalized two loop soft subprocess s^^^(A^o<3/(^A*);"s(/^)) 

(3.55) 



where we find 

^2- 



- 9 



^(2) 



(3.56) 



^ ^ /I 5 28 TT^ 

11 , 1/67 Tr^\ , 11 2 . 101 7 



5 777r2 ^ 13C3 



81 216 18 , 

/II 1 / 67 7r2\ 1/101 11 2 7 
+ ^^^^Vl6^"?(,144"48j"7Vl08" 288'' " 8^^ 

871.^+!![! + ll^y (3.53) 



324 864 120 72 



Renormalization properties of Wilson loops have been studied in detail in [69-71]. The 
Wilson path we considered has two cusps and light-cone segments leading to additional light 
cone singularities. This leads us to the following evolution equation [69, 70] 

(^/.|- + /?(a,)^) log 5 (^^,«,(^)j = -2rcusp(a.(/x))log^ - 2r,oft(a.(/x)), 

(3.59) 

where rcusp(5') is the well-known universal cusp anomalous dimension while rsoft(fi') is a 
path-dependent coefficient often called soft anomalous dimension. The factor 2 in front of 
the cusp anomalous dimension in the evolution equation (3.59) counts the number of cusps 
in the integration path. 
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The two quantities rcusp(5') and Fsoftly) can be evaluated through 0{a1) considering 
the 0{as) and 0{af) counter-terms 



ii'/=C;,i^^ + ilog^), (3.60) 



1 5 1/ 7 TT^ 





^ 7 

V 27 


' 67 




144 


48 y 



4e3 36e2 eV 27 72^ 

if 

e V 108 288 



^ ^ f 11 1/67 iT^\ 1 /lOl 11 2 ■ > > > > 



leading to the following results 



hO(a2). 

(3.63) 

The two-loop cusp anomalous dimension (3.62) was computed in [72], while the two-loop 
value of rsoft(9') was first deduced in [42] using renormalization group invariance of the 
cross section but it was never obtained by a direct calculation. 

Exploiting the non-abelian exponentiation theorem (3.11) and (3.12), we derive a com- 
plete two-loop expression for the soft subprocess 

.(^,..,.,).i-c.^(^..V||) 



I -^^2 ^ ^ /67 7r2\ 10^ ^ \ 2 ^/^ 
+ [^C^F- CaCf[- --]+ -CfTfuf) log^ ^ 



28 7r2\ /II o 101 7.\\, Nii\ 



C.T,„, - _ - _J + C.C. (^-.^ + - - -6J J log ^ 

+ Sf + 0(af), (3.64) 

where ^q^^ is the non-logarithmic piece at two-loop order for which we provide an analytic 
expression 

777r^ 13C3 



216 18 



^ ^ V 324 864 120 72 y y ^ ' 

The constant part (3.65) has also been calculated as a specific case of the two-loop soft 
hemisphere function in a work done in parallel with ours [52], in full agreement with our 
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result. Previously, it had been fitted by two different groups [42,44,51] using the Monte 
Carlo program EVENT2 [15]. Their results are reported below 

^(2), [42] ^ ^1^((58 ± 2)C| - (60 ± 1)CfCa + (43 ± l)CFTFnF), (3.66) 

^(2), [44] ^ ^00 (48.7045CI - {57.8)CfCa + {A3A)CFTFnF) , (3.67) 

S^^^' = (48.7045C| - (58.8 ± 2.25)CfCa + (43.8 ± 3.06)CFTFnF) , (3.68) 

while we obtain 

= (48.7045C|, - (56.4989)CaCf + (43.3905)CFrFnF) , (3.69) 

which is partly consistent with [44,51] but not with the earlier superseded numbers of [42], 

~(2) 

with the exception of the CfTfUf term. Notice that the determination of the constant 
is relevant to the matching of NNLL resummed cross section to the NNLO one since it is 
part of the G31 coefficient as it will be shown in the next section. 

We now solve the evolution equation for the soft subprocess (3.59), yielding 

siZsi^)) = " ^ (2rcusp(a.(fe^)) log ^ + 2r,Mk'))) , (3.70) 
where is the renormalization scale of the process. We can perform the substitution 

and exchange the order of integration in the integral containing Fcusp getting 

(3.72) 

We set fiR = Q in order to minimize the logarithmic corrections coming from the hard 
function in (2.28), moreover we replace fi^/Q^ = u'^ in the first integral of Eq. (3.72) and 
k^/Q^ = in the second one and we finally obtain 



(3.73) 



5(1, a,(^)) ^ Ju^Q 



The thrust observable is symmetrical under the exchange of the two hemispheres, so 
the factor 2 in the exponent accounts for the identical contributions due to both of them. 
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4. Resummation of large logarithms 



In the present section we derive a resummed expression for the cross section (2.13) in 
Laplace space starting from the renormahzation group (RG) evolution of each of the sub- 
processes. The effect of soft gluons has been taken into account in the previous section, 
but we still need to consider logarithmically enhanced terms due to hard gluons moving 
collinearly to one of the hard quark legs. Such an effect is encoded in the jet subprocess 
which describes the decay of a hard quark into a jet of collinear particles. The same sub- 
process can be found in other relevant QCD processes such as deep inelastic scattering and 
5-meson decay and it obeys the following evolution equation [49,58]: 

('^^ "^^ ^^^"^ (VW^'""*^^-*) " 2rcusp(a5(/^))log-^^ - 2reoii(Q;s(M))- 

(4.1) 

The collinear subprocess can be defined as a cut propagator of a massless quark in the 
axial gauge [37]. Indeed, the factorization used here is manifest in the axial gauge [37]. 

Equation (4.1) can be solved following the same technique used with the soft subprocess 
but now replacing fJ.'^/Q'^ = u and k"^ /Q"^ = u in the first and second integral in the exponent 
respectively. It leads to 



j(l,a.(ytQ)) 

(4.2) 

We now combine (3.73) and (4.2) together in the expression of the cross section 
^n{Q\ cXs) = H{1, a.ifx = Q))P asiQ)) S f ^, a,(Q)') = 



AT 



2 fl — (rsoft(a.(u2Q')) + T^oixMuQ^)))] , (4.3) 



where we explicitly set the renormahzation scale fj, = Q. 
Using the relation 



r^oitMu^Q')) = T^^UMuQ')) - ^Piasie) f^^°'f'^'''^\ (4.4) 

^u2Q2 k^ das 
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we recast Eq. (4.3) as 



aN{Q'',as) =H{l,as{Q))J' ( l,a,(A/^Q) 1 S [ l,a,(-^°^ 




2 



-'k^i L' ir-^^<°'<*='» |, (4.5) 



where we defined 



Arias) = rcusp(«.) - ^(«.)^%^, 

Br (as) = rsoft(as) + rcoii(as)- (4.6) 

The two coefficients Arids) and Br (as) can be computed in perturbative QCD. To this 
end we observe that the AltarelH-Parisi splitting function Pqq(as,z) fuffils the following 
limit [73] as z — )■ 1 

Pqqias,z) = 2^^!!^ + 2Br{as)S{l - z) + (4.7) 

where the dots stand for regular terms in the z — t- 1 limit. The asymptotic expression (4.7) 
is valid to all orders in perturbative QCD and it can be easily proven in the context of 
deep inelastic scattering as shown in [50]. The Mellin transform of the structure function 
F2{Q'^,x) can be indeed factorized in the threshold limit x — > 1 as a product of a hard 
virtual function H{Q/iJ,,as{fJ')), a coUinear jet function and a parton distribution function 
(f)q{N,^). Both the hard and collinear jet functions arc essentially the same ones as in the 
thrust case (up to constants in the hard subproccss due to crossing). The collinear jet 
function evolution is described by Eq. (4.1), while the hard function RG equation reads 



2 

+ 2(rsoft(a.(/x)) + 2r,oii{as{n))), (4.8) 



and the parton distribution function evolves according to the Altarelli-Parisi equation [74] 
(^^1^ + ^(a,)^^ l0g(l>q{N,fi) = Pqq{as,N), (4.9) 

where 

Pqqias,N) = - [ dzZ^-^Pqqias,z) (4.10) 
^0 

is the Mellin transform of the splitting function using the conventions of [75,76]. We now 
observe that for the structure function F2{Q'^,x) to be RG invariant, we have to require 
that the anomalous dimensions of the hard, jet and parton distribution functions sum up 
to zero, proving Eq. (4.7). 
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As stated at the end of the previous section, since the thrust is symmetrical under the 
exchange of the two hemispheres, we can factorize the soft subprocess as a product of two in- 
dependent "hemisphere" soft subprocesses. It follows that the resummed cross section (4.5) 
can be recast as a constant term multiplied by the evolution of two independent jets each 
of which is the product of the collinear jet and the respective "hemisphere" soft subprocess. 
This is in analogy to the structure obtained at NLL using the coherent branching algo- 
rithm [16]. We compare the two expressions observing that the only difference between the 
two exponents is the term l3{as)dTi,oft{ois)/das which gives a non- vanishing contribution 
only beyond NLL. It essentially accounts for large angle soft emissions whose effects do not 
contribute at NLL (it is easy to see that the first non-trivial term arises at 0{a^)). A second 
interesting feature which shows up beyond the NLL approximation is the interplay between 
constant terms and logarithms due to the factor J^(l, as{^/ {NQ/N)Q))S{l,as{NoQ/N)) 
that will be analyzed below. 

The resummed cross section in the dijet limit (4.5) takes the following form 

^iv(Q',a.)= + )sA^(Q^a.), 

SA^(Q^a.) = exp|L/i (^/3oL) + (^/3oL) + ^/3o/3 (^/3oL) + G31 (g)' L 



+ 0(«^l2)|, (4.11) 

where here L = logN. The function fi{(as/TT)f3oL) resums all the leading logarithmic 
contributions q^L"^"'^, f2{{as/TT)f3oL) resums the next to leading terms a^L" and so on. 
We furthermore require that fi{0) = so that at N"LL we can write 

k 



(4.12) 



fe>i 

k 



fn+1 {^PoL) = GuMi-n (^) ri>\. (4.13) 

^ k>n ^ ^ 

With this notation we see that the term GsiafL is a N^LL contribution due to the Taylor 
expansion of /4(^/3o-L). Nevertheless, such a term is relevant for the i?-matching of the 
NNLL resummed cross section to the ©(a^) fixed-order result, which will be discussed 
below. After expanding the functions At{ois) and By^ols) as 

Arias) = Y: A^'^ {-) \ Br{as) = ^ B^'^ (^) (4.14) 

in Eq. (4.5), we observe that A^^^ gives rise to terms of order q;"L"+^~'^ while B^^^ con- 
tributes with terms of order a^L^^^"*^ with n > k. 

The previous property ensures that the knowledge of B^^^ is sufficient to compute G^i 
and we do not need to know A^^^i which has not been computed yet (Tsoft is known at three 
loops [42], but the four loop value of Tcusp is still unknown). 
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The coefficient H{1, as{Q)),P{l, as{^/iNo/N)Q))S{l, as{NoQ/N)) in Eq. (4.5) con- 
tains all the constant terms and it can be expanded in perturbation theory. The function 
H{l,as{Q)) is known at three loop order [47,48] (this result has subsequently to be nor- 
malized to the total hadronic cross section a) and the two-loop non-logarithmic value of 
the collinear subprocess was computed in [49]. The constant part of the two-loop soft 
subprocess was evaluated in the previous section. We see that the coupling is evaluated 
at different scales in each of the three functions. We use the expression for the running 
coupling (A. 3) to express them in terms of as{Q) evaluated at the renormalization scale 
fJ-R = Q- The resulting expression has additional resummed logarithms sitting outside 
the exponent of Eq. (4.5) due to the running of and giving a well defined and finite 
contribution at large N. Such terms contribute from NNLL on and do not exponentiate 
naturally. Nevertheless, in order to bring the cross section to the form (4.11), we raise 
them to the exponent and we expand them to the desired order. 

One finds that the one-loop constants of the collinear and soft subprocesses contribute 
to f3{{as/7r)f3oL) while the two-loop ones contribute to G31. 

To evaluate the integrals in Eq. (4.5) we use the renormalization group equation (A.l) to 
change the integration variable to ag. After imposing the normalization condition /i(0) = 
we find 



/i(A) 



4^[(1 - 2A) log(l - 2A) - 2(1 - A) log(l - A)], 



(4.15) 



'0 



/2(A) 




^(2) 5(1) 

-^[21og(l - A) - log(l - 2A)] + 2-- log(l - A) 



^[log(l - 2A) + ^ log2(l - 2A) - log(l - A)(2 + log(l 



A)] 



(4.16) 



/3(A) 



2c: 



■'^^ A 2cf A 25(2) ^ ^(3) ^ 



/3o 1-2A" /3o 1-A I3i 1-X (1 - A)(l - 2A) 

2A^^^-fE A A^'^^Pi 3A2 + (1 - A) log(l - 2A) - 2(1 - 2A) log(l - A) 

/32 (1-A)(1-2A)+^ (1-A)(1-2A) 

a 25«/3i a + log(l - A) 




+ 



(1 - 2A) log(l - A)] + g[-A2 + (1 - 3A + 2A2)(21og(l - A) - log(l - 2A))]] 





[A2-(l-2A)log(l-A)(2A + log(l-A))] , (4.17) 
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where A = ^^/3ologiV. 

At NLL the previous functions reproduce the result obtained in [16]. The constants cj^ 
and are the one-loop non-logarithmic terms of the collincar and the soft subproccsses 
respectively. They arise from the term j2(l, a,(^(iVo/iV)Q))5(l, as(iVog/iV)), as dis- 
cussed above. The coefficients A^'^ and B^'^ can be determined as shown in Eq. (4.6) using 
the two-loop value of VsoU computed above and the three-loop splitting functions [75,76]. 
They are reported in Appendix A. 

We observe that the normalization condition /t(0) = is automatically fulfilled by 
both /i(A) and /2(A) while it has to be imposed to obtain /3(A) (4.17). This could be 
considered as a signal of the breakdown of natural exponentiation beyond NLL. Forcing 
such a constraint gives rise to a residual constant value which has to be taken out of the 
exponent and that contributes to the constants We will determine the value of such 
constants directly in thrust space in section 4.1. 

4.1 Inversion of the integral transform 

In the present section we perform the inverse Laplace transform of the resummed cross 
section. We recall the definition of the normalized cross section 

Rt{t) = 2^ ^e^"c^A^(Q^ as) + 0(r), (4.18) 

where criv(<3^, c^s) is defined in Eq. (4.11). The contour C runs parallel to the imaginary 
axis to the right of all the singularities of the integrand function. A method to invert (4.18) 
was proposed in [16]. Keeping their notation we rewrite Sjv(Q^,as) in (4.11) as 

Siv(Q^as) = e-^("^(«')'^°sAf)^ (4.19) 

where 

P{as{Q\L) =Lh ^oL) +/2 {^^oL) + ^Pof3 {^|3oL)+G^^(^^yL 

+ 0{aiL'^). (4.20) 
With this notation CTiv((5^,«s) reads 

aN{Q',as) = (l + f;C.(^)')e^(-(«^)>^^). (4.21) 

We now Taylor expand the exponent T with respect to log A?^ around log A?^ = log(l/r) up 
to the desired order and we end up with the following expression for Rt{t) 

Rrir) = (1 + E (i) ') e^^-^'^^)'^"^ ^ |^ ^exp [z. + T^') {as{Q% log 1) log ^ 

(4.22) 



+ (a,(Q2), log i) log^ ^ + 0« log"-2(l/r)) 



-28- 



where v = Nt and 



Jf(-) (as{Q^),log- 



-F{as{Q')^og{l/T)). 



(4.23) 



TJ 91og"i 
More precisely, we have 

^«(«.(Q2),iog(l/r)) = /i(A) + A/((A) + ^/3o/^(A) + 0(a^log'^-2(l/r)), (4.24) 

vr 

(a,(Q2)jog(l/T)) = 2^/?o/((A) + ^/3oA/r(A) + 0«log"-2(l/r)), 



TT 



TT 



(4.25) 



where, from now on, A = (Q;5/7r)/3o log(l/r). Without loss of generality we expand the 
factor 



exp 



1 ^ 



J-(2)(a,(Q2)^log(l/^))log2^ 



(4.26) 



in a Taylor series in ^^^^ itself and we see that the series of terms [T^"^^)^ with k>2 gives 
at most rise to terms of order ©(a^ (ttg log(l/T))"), so we drop them. Using the result [16] 



1 f diy 

- — / — log*^ u exp v + ^^^^ log v 
2m Jc I' L 

we obtain 

i^T(r) = f 1 + f; C Jf^) \,-^(- W^)>^'>^ ^) 



^(^(i))fer(i-^{i))' 



(4.27) 



k=l 



^2-^ d72(A)r(l-7(A)) 



Lr(l-7(A)) TT 



+ 0(af(a,log(l/r))") 



+ ^M^(A)^^°^(^--(^)^ 



r(i-7(A)) 

(4.28) 



where wc neglected subleading terms and we defined 7(A) = /i(A) + A/((A). t/j^""^ is the 
(m + l)-th derivative of the logarithm of the Euler F-function. 

After replacing the coefficients and G31 in the previous expression with Ck and G31 
computed as explained in Appendix B, we can recast Rt{t) as 



i2r(r)=(l + f]Cfc 



k=l 



as 
27r 



3.^(a.(Q^),logl) 1 1 + ^^o/^(A)^(°)(l - 7(A)) 

i(l-7(A)j[ TT 



a,. 



+ — I 7£ I 2 - 7e J - y 



TT 



(4.29) 



where now 

J^(«,(Q2),l) =L/i ('^fi^L)+h (-PoL) + ^^0/3 (^^oi^) +G31 
+ 0(afL2). 



^"l L 
2tx J 



(4.30) 



The expression (4.29) for Rt{t) resums all large logarithms through NNLL and it holds 
up to terms of order ©(a^ log^(l/r)). 
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The constant term agjiiCp (7£;(3/2 — ^e) — ^^/6) is defined such that the term in 
square brackets is normahzed to one as A — > and ensures that all the constant terms are 
included in the Cfe coefficients which read 



C\ =Cf 



5 TT^ 

2 + T 



41 



C2 =CU^- o + ^TT* - 6C3 + CaCf 



36 



/ 491 



537r2 IItt'' 104C3 



108 



360 



28C3 



6 27 

The coefficient G31 is found to be 



G: 



31 



77 



36 ^ 45 



Stt^ 104C3V^^ / 



191 eivr^ 



+ cl 



81 

leiTT^ 



+ 



72 



32C3 

9 

537r^ 



+ 



9 

+ C\Cf 
452C3 



1118 6447r^ 



81 



17V; 
135 



(4.31) 



(4.32) 



292C3 



/ 5951 



66257r2 3837r^ 404(3 

+ -^TTT- + + IOC5 



324 



540 



9 



9 



+ 27r"C3 + 3OC5 



29 57r2 



+ 53(3 - ^T^Cs + 132C5 
15 3 



(4.33) 



The value of the 0{a^g) constant C3 is unknown and we fit it using the Monte Carlo 
program EERAD3 [21]. Details about the fitting procedure are explained in the next section. 
Following the conventions of [16] the final resummed expression Rt{t) (4.29) can be written 
as 



RT{T)=[l + Y,Ck 



k=l 



a 

2lT J 



1 lexp 



log -gi{\) + 52 (A) + ^/3o53(A) + ( ^ ] G31 log 



/a£\3^ . 1 



(4.34) 



where 



ffi(A)= /i(A), 

ff2(A) = /2(A) - logr(l - /i(A) - A/{(A)), 

53(A) = /3(A) + [j[ + ^A/f (A)) - 7(A))' - ^«(1 - 7(A))) + /^(A)V'(°)(1 - 7(A)) 

Cf 



+ 



7£ o - 7E - 



TT 



(4.35) 



notice that the functions gi(X) do not generate any constants (i.e. ^1(0) = 0). 

The resummed expression at different logarithmic orders evaluated around the peak 
region is shown in Figure 6. 

Wc observe that an exact inversion of Eq. (4.18) requires to choose the integration 
contour C to the right of the Landau singularities present in the resummed functions fi{\). 
Following our prescription [16] for the inversion, we avoid such singularities by expanding 
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Figure 6: Comparison of the resummed result at different logarithmic orders around the peak 
region. 

around log N = log ^ and then integrating by means of the residue theorem closing the 
contour in the left half-plane using Eq. (4.27). In doing so, we are neglecting the contribu- 
tion due to the residue at the Landau pole, which gives rise to power suppressed terms [77]. 
The Landau singularity is then mapped onto the thrust space without contributing to the 
Laplace inversion. By expressing the SCET result of [42] in the form (4.34) [78], we obtain 
full analytic agreement on (73(A). This is a non-trivial result, since in [42] the scales in 
the resummation kernels are fixed at the outset of the calculation, while in our analysis, 
they are integrated over in the Laplace inversion. This difference results in a different 
treatment of the Landau singularity in the resummed expressions, which could produce a 
power-suppressed difference between the cross sections obtained in both approaches. The 
exact agreement of our 53 (A) with [42] demonstrates that there is no Landau-pole ambigu- 
ity between the two approaches. However, in [42] subleading terms arising from the scale 
fixing in T-space and part of the constants are kept in the exponent resulting in a numerical 
difference when compared to our result. 

4.2 Determination of the 0(0^) constant C3 

In order to match the resummed result to the NNLO cross section the 0{a^) constant C3 
must be extracted from fixed-order data. We do it by subtracting the logarithms from 
the fixed-order C(t) coefficient obtained from EERAD3. The logarithmic part ^iog(''") of the 
NNLO coefficient is decomposed into its different color contributions according to 



: r{3) 



, 7?{3) 
^0 + ^log 



1/Af2 



+ R 



np N 



+ R 



np/N 



+ R 



(3) 



(4.36) 



EERAD3 is run with a technical cutoff yo = 10~^ which affects the thrust distribution below 
To ~ \/yo- This forbids us from probing the far infrared region and we perform the fit for 
values of r larger than tq. Numerical fixed order results are obtained with 6 x 10'' points for 
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Figure 7: Fit of tlie 0{a^) coefficient C3. The blue band shows the statistical error of C3 and the 
red arrow indicates the fitting interval. 

the leading colour contribution and 10^ points for the subleading colour structures. Because 
of the presence of large fluctuations in the Monte Carlo results, each color contribution is 
fitted separately over an interval where the distribution is stable and the different results 
are combined to find the numerical value of C3. The results of the fits and the different 
fit intervals are given in Table 1. As an alternative approach we first sum all the color 



Color 


(-logr)min 


(-logr)max 


Fit result 




4.2 


5.2 


3541 ± 51 




4.2 


5.4 


-265 ± 8 




3.8 


5.2 


-71 ±3 


Nnp 


4.6 


5.6 


-5078 ± 145 


np/N 


4.6 


5.8 


236 ± 7 


Up 


4.2 


5.2 


95 ± 120 


Sum of all colors 


-1543 ± 147 


All colors 


4.2 


5.2 


-1051 ± 178 



Table 1: Intervals and results of the fits for C3 for the different color contributions. 

contributions to the C(r) coefficient, then we subtract Eq. (2.17) and finally fit C3. The 
result of the second approach is shown in Figure 7. We consider the difference between the 
two approaches as a systematic error and as final result we obtain 

C3 = -1050 ± 180(stat.) ± 500(syst.) . (4.37) 

Considering that there is no statistical correlation between different bin errors, as a different 
possible estimate of the systematic uncertainty due to the sizeable fluctuation, we varied 
the fit range observing that it does not alter the result in any significant way outside the 
quoted systematic error margins. It is worth stressing that the numerical impact of C3 on 
the distributions is negligible, as it will be shown in Section 6, such that the large relative 
error range is tolerable for all practical purposes. With the determination of the 0{a^) 



-32- 



constant we now have all the needed ingredients to perform the matching of the NNLL 
resummed result to the fixed NNLO. 

5. Matching of resummation to fixed-order calculations 

There are different matching schemes proposed in the literature, however mainly two are 
used: the i?-matching scheme and the log(i?)-matching scheme [16,19]. The new results 
presented in the previous sections allow us to compare for the first time the predictions of 
the two schemes at NNLL+NNLO accuracy. In the i?-matching scheme the two expres- 
sions (2.11) and (2.13) are matched and logarithms appearing twice are subtracted. The 
explicit expression for the matched integrated cross section R{t, Q) depends on both the 
logarithmic and fixed-order accuracy considered in the matching. At NNLL+NNLO the 
following formula holds (for the sake of clarity wc drop any dependence on the renormal- 
ization scale, which will be analyzed separately and write only UgL as arguments of the gi 
functions) : 

R{t,Q) = (1 + Cias + C2al + Gso^) e(i9i(°=i)+92(°=i)+^/3o93(a.L)+«!LG3i) 

+ as (Air) - (r)) + (b{t) - itig (r)) + (c(r) - itlg (r)) . 

(5.1) 

The terms in the second line correspond to the remainder functions (ij(r) defined above. 
It is however preferable to write it as difference between the full fixed-order coefficient and 
its logarithmic part since these are the functions which are known in practice. 

The log(i?)-matching scheme [16] is believed to be theoretically the most stable one 
and for this reason it is generally preferred [19]. In this case the matching procedure is 
given by 

log (i2(r, as)) =Lgi{asL) + g2{asL) + aggsiagL) 
+ as {A{t)-GuL-Gi2L^) + 

+ al (b (r) - ^Al (r) - G21L - - ^23^^=^) (5-2) 

+ al (c (r) -A{t)B (t) + ^A^ (r) - G32L'' - G^^L^ - G^L^ 

It is worth noting that the dependence on the 0{a^g) coefficients C3 and G31 disap- 
pears in this scheme. The matching procedures presented above are valid over the whole 
phase space. However, unlike to the fixed-order prediction (2.11), in which every coefficient 
vanishes in the kinematical limit such that i2(rmax) = 1, in the limit r Tmax the two 
predictions (5.1) and (5.2) give a wrong {i.e. non-vanishing) result. A fixed-order calcula- 
tion takes into account only a finite (in fact very few) number of final state particles so the 
differential cross section d^j j (\t obviously has to vanish at the kinematical limit t — '7~max* 

For the NNLO fixed-order prediction, the maximum number of final state jets is five. 
Therefore the cross section should vanish at the kinematical limit for six partons. The 
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limited predictive ability of the two matching schemes in the multijet region can be solved 
by modifying them slightly. This is done by imposing a kinematical constraint, which 
assures the right prediction for r Tmax- 

The constraints for the so-called modified log(i?)-matching scheme are [19] 



log (E(Tinax, as)) = 



1 da(r) 
a dr 



di?(r) 



dr 



0. 



(5.3) 



In order to fulfil these two constraints, we follow the prescription proposed in [19] and 
redefine L as follows: 

IV / 1 ^^ 



L 



1 
P 



log 



+ 1 



(5.4) 



The power p is called "degree of modification". We choose p = 1, as usual in literature. 
It determines how fast the integrated cross section is damped at the kinematical limit. 
The value of rmax is given by symmetry arguments and at LO and NLO can be computed 
exactly giving respectively rmax.LO = 1/3 and rmax.NLO = 1 ~ l/\/3- At NNLO we can fix 
it using the result given by EERAD3 at r^ax = 

0.4275. 

For the log(i?)-matching scheme, the substitution (5.4) is sufficient to fulfil the con- 
straints (5.3). In the i2-matching scheme one further modification is needed 

p- 

1-1 — I , i = 1,2,3, (5.5) 



Ga^Gair) = Ga 1- ( , ^ = 1,2,3, 

leading to the following expression for the /^-matching at NNLL+NNLO: 
R (r) = (1 + Cias + C2a1 + Cga^) 



X exp 



L'giiasL') + g2{asL') + ^Pogsia^L') + a^L' G3i(r) - 

TT 

+ as [A{t) - b!^^ (r) |L^L',Gii^Gii(r)) 

+ «s - ^log i'^) k-^L',Gii->Gii{T),G2i->G2i(T)) 

+ al (c{t) - b!S(t)\l^l' 



(Gil + G21) 



7,Gii^Gii(t),G21^G2i(t),G31^G3i(t) 



(5.6) 



The dependence on the renormalization scale was so far not considered. Every term 
beyond the leading order acquires an explicit //-dependence, which for the fixed-order 
coefficients is given in (2.4). For the resummation functions gi (agL) the renormalization 
scale dependence is given by 



g2{asL) g2 iagL,!?) =g2{asL) -\ {asLfg[{asL) logixl) 



gsiagL) gs {asL,^x^) =g3{asL) + 

+ 



{asL)g2{asL) + -^{asLfg'i{asL) 

■n-po 



^{a,Lfg[{a,L) + ^{asLfg'HasL) 

TT ZTT 



(5.7) 

log(a;^) 
log(4)'' 



(5.8) 
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where g[{asL) stands for the derivative of gi{asL) with respect to a^L. Correspondingly 
the coefficients Gij and Ci change as foUows: 



G21 - 
G22 - 
G33 - 

G31 - 
G32 - 
C2 
C3 



G21 
G22 
G33 
G31 
G32 

y C2 
> C3 



G21 + 2/3oGii log(4) 
G22 + 2/3oGi2log(x2) 



G 



33 + 2/3o2G23log(x2 



G31 + (2/3olog(x2))'Gii + 21og(x2)(2/3oG2i + 2/3iGii 

G32 + (2/3olog(4))'Gi2 + 21og(4)(2/3oG22 + 2/3iGi2 
G2 + 2/3oGilog(4) 

G3 + 2 (2;3oG2 + 2AGi) log(4) + (2^0 log(4))' Gi. 



(5.9) 



One further source of arbitrariness is the choice of the logarithm to be resummed. 
In fact, it is not clear whether powers of Q;slog(l/T) or powers of aslog(2/T) have to be 

resummed. The origin of this arbitrariness has to do with how much of the non-logarithmic 
part of the fixed-order prediction is exponentiated together with the logarithms. We can 
express this arbitrariness by introducing a new constant xl, which rescales the logarithm 
to be resummed [19]: 



L ^ L 



1 
P 



lyXLTj yXLTmaxJ 



+ 1 



(5.10) 



This rescaling modifies once more the resummed formulae and their expansion coefficients. 
By requiring R{t) = R{t) , where i?(r) denotes the rescaled integrated cross section 
according to (5.10), we find the following replacements: 



Gi(r) = Gi + Gn log(xL) + G12 log{xLf , 

G2(r) = G2 + (G21 + GiGii) \og{xL) + (g22 + ^Gj, + G1G12 ) log{xLf 

+ (G23 + G12G11) \0g{xLf + ^Gl2 log(xL)^ , 

G3(r) = G3 + (G31 + G1G21 + G2G11) log(xL) 

+ (g32 + G1G22 + ^GiGfi + G2G12 + GiiG21^ logixLf 

+ (g33 + G11G22 + G12G21 + G1G11G12 + ^G?i + G1G23) log(xL)' 

+ (g3a + G12G22 + lciGl2 + G11G23 + ^GfiGia) log(xL)^ 
+ (G12G23 + IgI^Gu) log(xi)^ + ^Gl^\og{xLf , 



(5.11) 



(5.12) 



(5.13) 
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Gl2 - 


G\2 


Gil - 


■> Gil 


G22 - 


G22 


G33 - 


■> G-iz 


G21 - 


■> G21 


Gz2 - 


■> G32 


G31 - 


■> G31 



G12 , G23 — G23 , G34 — G34 , 
Gil + 2Gi2 log(xL) , 

G22 + 3G23 log(xL) , 
G33 + 4G34 log(xL), 

G21 + 2G22 log(xL) + 3G23 log(xL)2 , 
G32 + 3G33 log(xL) + 6G34 log(xL)^ , 

G31 + 2G32 log(xL) + 3G33 \ozixLf + 4G34 \oz{xLf . (5.14) 
The corresponding changes of the gi functions axe 

g\{asV) g\{asL) =gi{asL), 

g2{asL) hicts 'L) =g2{asL) + [gi{asL) + {as 'L)g[{as 'L)^ log(xL) , 

9z{oisL) gsioisL) =g3{asL) + —g2{asL)\og{xL) 

Po 

+ (^g[{asL) + ^iasL)g';iasL)^ logixL^ • (5.15) 

The transformations due to a variation of and xl are completely general and hold for 
all possible event-shape obscrvablcs which can be described with this matching formalism. 
Furthermore the order in which they are carried out is not important since they commute. 

6. Results 

Having set up the matching formalism in a way to access the theoretical uncertainties, we 
can apply it to the case of thrust derived in Section 4 using the fixed order result from 
EERAD3 [22,23]. 

In Figure 8 we compare the weighted cross section of the new matched NNLL+NNLO 
results with the old NLL+NNLO derived in [32] . The modification due to the resummation 

is sizable, leading to a 8% increase of the distribution around the peak region. The cff'cct 
of the additional resummed subleading logarithms becomes progressively less important 
towards the multijet region, where the increase is nevertheless of about 5%. It is interesting 
to note that the matching of NNLO with NNLL resummation shifts the pure NNLO result 
also in the multijet region (Figure 9). This was not the case for NLL+NNLO, for which 
the impact of resummation in the region of large r was negligible. This is another sign of 
the importance of the NNLL contribution. 

The renormalization scale dependence, which was observed to increase from pure 
NNLO to NLL+NNLO [32,35] because of a mismatch in the cancellation of renormalization 
scale logarithms, is obtained by varying 0.5 < < 2. It decreases at NNLL+NNLO by 
20% in the peak region compared to NLL+NNLO. The magnitude of the scale uncertainty 
varies between 4% in the 3-jet region and 5% around the peak. 

In Figure 10 we compare the unweighted parton-level cross section at NNLO with the 
matched NLL+NNLO and NNLL+NNLO cross sections and with experimental hadron- 
level data from the ALEPH experiment [1]. We note that there is a visible shift of the 
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0.1 0.2 0.3 0.4 

1-T 



Figure 8: Comparison of the weighted cross section in the log(i?)-matching using NNLL+NNLO 
and NLL+NNLO. The plot on the top shows the two distributions, with the uncertainty band due 
to scale dependence. The curve in the middle shows the difference between NNLL+NNLO and 
NLL+NNLO normalized to the NLL+NNLO curve. The impact of the resummation at NNLL is 
an increase in the distribution of order 5-8%. The lowest plot shows the absolute scale dependence 
of the two curves. 

theoretical prediction towards the experimental data, which are best described by the newly 
computed NNLL+NNLO distribution. Around the peak region, where non-perturbative 
hadronization corrections are large, the parton level prediction fails to describe the data. 
Hadronization effects can account for this discrepancy. We will address this issue in a 
future publication. 

The computation of the two-loop constant C2 presented in Section 3 and the fit of the 
©(Og) constant C3 of Section 4.2 allow us to perform for the first time the matching in the 
i?-scheme using NNLL accuracy and fixed NNLO results. On the left plot in Figure 11 we 
compare the i?-matching and the log(-R)-matching scheme predictions at NNLL+NNLO. 
The difference between the two matching prescriptions is very small and lies well below 
the scale uncertainty. This implies a very good stability of the theoretical predictions 
under variation of the matching scheme. Because of the big uncertainty in the value of the 
constant C3, we vary it within its error to investigate the phenomenological impact on the 
cross section. The results are shown in the right plot of Figure 11. The upper plot shows 
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Figure 9: Comparison of the weighted 
cross section at NNLO with the matched 
NNLL+NNLO predictions. The contribution 
of NNLL resummation is sizable over the full 
thrust range. 



NNLL+NNLO 
NLL+NNLO 
NNLO 




Figure 10: Comparison of the parton 
level fixed-order prediction at NNLO to 
NLL+NNLO and NNLL+NNLO distributions. 
Experimental data are taken from the ALEPH 
collaboration [1]. 
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Figure 11: The left plot shows the comparison of the i?- matching scheme and the log(i?)-niatching 
scheme results. The width of the curve shows the uncertainty related to the scale variation. The two 
matching schemes agree very well over the full thrust range. The right plot shows the impact of the 
variation of the 0{al) constant C3 on the distribution in the R-matching scheme. The difference 
is at the per mille level. 



the different distributions obtained by setting C3 = —500, —1000, —1500 respectively. The 
three curves are almost indistinguishable and the tiny fluctuations in the distributions are 
due to the fluctuations of the NNLO result. In the lower plot we take the distribution with 
C3 = —1000 as reference and plot the difference between the reference cross section and 
those obtained with C3 = —500 and C3 = —1500 respectively. The difference is less than 
1.5%o and it is therefore completely negligible compared the other theoretical uncertainties. 

As discussed in the previous section, another source of theoretical uncertainty is the 
choice of the logarithms to be resummed. We can estimate this uncertainty by varying the 
parameter xl- In Ref. [19] several prescriptions are given on how to set the correct variation 
range for xl for different observables. For the sake of simplicity and since we are not 
performing a fit of the strong coupling constant, we choose to vary xl within the canonical 
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R-matching NNLL+NNLO ■ 
Log(R)- matching NNLL+NNLO ■ 




Figure 12: Dependence on the resummed logarithms, determined by varying the parameter x^. 
The left plot shows the change in the xl dependence between NLL+NNLO and NNLL+NNLO. 
The upper plot shows the distributions with the corresponding uncertainty band, in the lower plot 
we compare only the uncertainties. In the right plots the x^ dependence using the two different 
matching schemes is shown. 



interval 0.5 < < 2, similarly to what is chosen to quantify the renormalization scale 
uncertainty. This choice is also close to the nominal range of variation proposed in [19]. The 
impact of this variation is shown in Figure 12. The left plots show a comparison of the xl- 
dependence between NLL+NNLO and NNLL+NNLO predictions. The lower plot allows to 
quantify the reduction of the uncertainty due to a variation xl . Apart from the far infrared 
region, it is observed to decrease by 50% in the peak region. The scale-dependence reduction 
is smaller towards the multijet region, where the contribution of the logarithmic part 
becomes less important. The resummation uncertainty at NNLL+NNLO varies between 
2% and 3%. In the right plots the same comparison is made at NNLL+NNLO using the 
i?-matching and log(ii)-matching schemes. We observe a similar x^^-dependence in both 
schemes. 

7. Conclusions and outlook 

In this paper, we computed the resummed thrust cross section in e~^e~ annihilation to 
NNLL accuracy in Laplace space and we provided a compact expression for it in momen- 
tum space. We analytically derived the two-loop corrections to the soft subprocess for 
thrust, for which only two numerical results [44,51] fitted to the subtracted fixed-order ex- 
pressions were available previously. We find reasonable agreement with them within their 
uncertainty. Our NNLL result for the thrust distribution confirms an earlier calculation 
carried out in the SCET framework [42] in a different calculational approach. We find 
full analytic agreement between the two results, including the absence of power-suppressed 
terms which would be beyond the logarithmic accuracy. 

The NNLL corrections were matched to the fixed-order NNLO prediction for the 
thrust distribution in the i2-matching and the log(i?)-matching schemes. By inspecting 
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renormalization-scalc logarithms in the resummcd and fixed-order expressions, the appro- 
priate combinations are NLL-^NLO and NNLL-^NNLO, while the NLL-hNNLO combina- 
tion leads to a mismatch in these logarithms, resulting in an enhanced scale-dependence. 

In SCET-based analyses [42-44] the appropriate combinations are usually determined 
by requiring that the order of ag of the hard, jet and soft functions matches the fixed order 
expansion. This results in a different identification of matched distributions. We observe 
that the NNLL contribution is sizable even outside the dijet region and that the inclusion 
of the NNLL corrections further stabilizes the perturbative prediction throughout the full 
kinematical range. 

A comparison with experimental data has to include hadronization effects, which are 
most pronounced in the two-jet region. In the past, these were often obtained using 
leading-logarithmic parton shower Monte Carlo programs, which turn out to be clearly 
insufficient [35] in view of the precision now attained by the perturbative description. Sys- 
tematic approaches to hadronization within the dispersive model [31, 79] or by using the 
shape function formalism [43] are offering a more reliable description. We will address this 
issue in future work. 



Acknowledgments 



We wish to thank L. Magnea for interesting discussions on factorization issues, L. Yang 
for stimulating comments on several aspects of the calculation and G. Heinrich for help 
concerning SecDec. Furthermore we acknowledge several fruitful discussions with A. Banff, 
T. Becher, E. Gardi, M. Grazzini and D. Maitre. GL would like to thank the University 
of Zurich for its hospitality. TG would like to thank the Kavli Institute for Theoreti- 
cal Physics (KITP) at UC Santa Barbara for hospitality while this work was completed. 
This research is supported in part by the Swiss National Science Foundation (SNF) un- 
der contract 200020-126691, by the UK STFC, by the European Commission through the 
"LHCPhenoNet" Initial Training Network PITN-GA-2010-264564 and by the National Sci- 
ence Foundation under grant NSF PHY05-51164. 



A. Constants and anomalous dimensions 



The renormalization group equation for the QCD coupling constant reads 

1^ (f^^ ft , (^) ^3 ....). (A.1, 
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where the coefficients of the P{as) functions are 

Pi =^Cl - ^CaTfUf - ^CfTfuf, 

325 2 5033 2857 
^^=3456^^- 1152"^ + 12^' 

1093 , / 50065 809 , \ n / 1078361 1627 \ 891, 149753 



186624 V41472 2592^ / V 41472 1728 V 64^ 256 

(A.2) 

Equation (A.l) can be solved in perturbation theory and it gives the following resummed 
expression for the strong coupling 

+ ^Itw (f ( ^ ^ - ^ - 1^^) 

- ^2A(l + A) + ^(2(1 + 2A)log(l + 2A) - 31og(l + 2A) + 2A(1 + 2A)))) 
Po Po / / 

+ ..., (A.3) 

where here A = {as{pR)/iT)l3Qlog{fi/ fin). The coefficients A^'''^ and i?*-*) used in the re- 
summed cross section (4.5) can be computed using Fsoft at two-loop order and the three- 
loop Altarelli-Parisi splitting functions computed in [75,76]. Using the notation (4.14) we 
find 

^« = Cf, 

^(3) ^ CFnl ] + CaCfuf \ + CluF h — 

^V324 2167 ^ ^ ^ \^ 1296 72 / ^ \ m 2) 

^^2^ / 15503 3897r2 llvr^ UCsX 
+ ^^^^ 1^^592" " ~86r + W " ' 

(A.4) 
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V 16 1296 2880 72 J ^\576 648 36 

2 /23 57r2 297r4 ITCsA ^ ^ ^2 / 151 2057r2 2477r4 2IIC3 



^64 288 4320 24 y ^ 256 576 8640 48 

1 2. 15C5\ ^2^ ^1657 28l7r2 ^4 97^^ 
" 24'' " ^ j + ^^^^ V2304 " "648~ + 1152 + ^ " X 
, r3 ^ 29 37r2 vr-i 17(3 .Is., 15(5 \ 



The functions H{l,as{Q)), J{l,as{\/{No/N)Q)) and S{l,as{NoQ/N)) can be expanded 
in a power series in the coupHng as 



Hil, am) = 1 + E4^^^, J (l> o.s{VWn)Q)) = 1 + ^ .(0 ^:(vWiV)Q) 
i>l ^ i>\ 



S{l,a.iNoQ/N)) = l + ^c«^^^^ll^. (A.6) 



1>1 



The coefficients cf^^ can be evaluated using the on-shell quark form factor [47,48], normal- 
ized to the total hadronic cross section 



4" -^4-? + ^). 



(2) ^2/745 137r2 677r4 15(3^^ ^/^5867 9l7r2 17C3\ 
4 -^H"6r""6- + l80--^J+^^"^^Hl296 -W-l^J 

V 5184 864 90 72 y ^ ^ 

The two loop non-logarithmic term of the collinear subprocess was computed in [49], re- 
sulting in 

II) „ /7 ir^ 



(2) 



/ 4057 137r2\ / 53129 1557r2 377r4 9(3 A 

^^^^"^ ^^Vl0368 " W" 2880 "Xj 

2^205 97^2 61vr4 SCsA 
+ ^^Vl28"T9^ + T440"Xy'' ^"^"^^ 

and the non-logarithmic part of the two loop soft subprocess was computed in section 2 
(3.64), (3.65) 



.(2) / 5 777r2 13C3\ ,^^f 535 871^^ 7^4 

- - 32^^ + ^^^^"^U + W-^J+^^^^V"^"^^ 



143C3 
72 
(A.9) 
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B. Computation of the G31 coefficient and of the constant terms 



The starting point to obtain Eq. (4.29) is the cross section in Laplace space (4.5) in which 
we include the coefficients vl^*-* and B^"^^ with i < 3. In Section 4 we show that this is 
sufficient to compute the N^LL function /3(A) and the N^LL coefficient G31. The latter 
combines with terms arising from the inverse Laplace transform to produce G31. With 
the normalization chosen for the functions /i(A) (Eqs. (4.15-4.17)) the constant terms in 
Eq. (4.11) are defined as 

Hi^,as{ii))exp (^^Gio + (^^20 + 0(a^)) =1 + ^C, + (^)' ^2 + ^(a^), 

(B.l) 

where Gio and G20 are the constant terms left in the exponent after computing /3(A) and 

^31. 

In the Laplace inversion we have to add the N^LL term g-^*-'^-* {as{Q'^), log ^) log'^ v in square 
brackets of Eq. (4.22) and include the next subleading term in the definition of J^^^^ and 
getting 

^W(a,(g2),iogi) = /i(A) + A/{(A) + ^Pom) + (^^0)' /3(A) + 0(a^ log""^ -), 
{as{Q% log -) = 2^M[{X) + ^/3oA/f (A) + (^/3o) ' /^(A) + 0« log"-^ -), 

T TT TT \ TT / T 

^^^\as{Q%log -) = 3 (^/3o)' /{'(A) + (^/3o)' A/r(A) + O^log^^ -). (B.2) 

It is important to notice that the terms /2(A) , /3(A) and the whole J^(^)(as(Q^),log ;^) 
contribute at most with logarithmic order C(a" log"~^ i) {i.e. N^LL) and we do not need 
them for a N^LL order resummation. Nevertheless, they are relevant for the computation 
of the coefficient G31, so we keep them. We now solve the integrals as shown in Section 
4.1 neglecting subleading (N^LL) terms and we expand the final result through order 
©(a^log^) obtaining Ci, C2 explicitly as reported in Eq. (4.31,4.32). To determine G31 
we observe that the coefficient of log ^ in the previous expansion is GS1 + C1G21 + C2G11, 
where G21 and Gu are known. The result is reported in Eq. (4.33). 

C. Evaluation of the integrals over soft gluons transverse momenta 

In the present section we show how to evaluate the integrals over the transverse component 
of soft gluon momenta used in the text. Let us consider the integral 



/ 



{q + fe)2 - iO 



Using the relations 



= {k- n){k ■n)-kl, d'^-^k± = _ ^^^ {kly'dkl, (C.2) 



7r2 



q^ = {q- n){q ■ n) - ql, d'^-^q^ = -^—-(ql)--dqlsm-^'{e)de, (C.3) 
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where 9 is the angle between the d — 2-dimensional euchdean vectors q± and k_\_ , we can 
recast (C.l) as follows 



TT^- vr5- /'^^2 ^^^(^2 • n){k ■ u) - kl)6M {{q ■ n){q ■ n) - ql) 



r(l-e)r(i-e)7 ^ v^^^ -2\k^\\qi_\ 

""cosio) - +io' 

(C.4) 

The integrals over q'^ and /c^ can be easily evaluated using the two 5 functions while the 
angular integral needs some attention. We set = jq ^^d we consider the 



angular part of Eq. (C.4) 



L 



cos(6') -JC + iQ 



d9, 



where it is straightforward to show that /C > 1. The previous integral can be evaluated 
by setting cos(^) — > 2i — 1 and using the integral representation of the Hypergeometric 
function 



cos(e)-/C + iO 1 + IC-iOJo 



T2(l-e) 2Fi(l,i-e,l-2e, 



2 



r(i-2e) i + x:-io • 



This leads to the solution of (C.4) 
J d^'-^q^d'^'^k^ 



,.-2„ ,.-2, <5(^W1H^ ^ VT^-- vrl- r^(^-6) 

{q + kf-m r(l-e)r(i-e) r(l - 2e) 



2Fi(^l,^-e,l-2e,4 
X {q ■ n)''^{q ■ fi)~^{k ■ n)'~^{k ■ n)~^ — 



(q-n){q-n)(k-n)(k-n) 



[^/{(f'n){k^ + ^/{q ■ n){k ■ n)) 

(C.6) 



used in the text. The second relevant integral is 

"^'^ ^^{{q+kr+miq+kr-ioy ^^-^^ 

appearing in the computation of the vacuum polarization diagrams. Such an integral can 
be evaluated easily using the previous reult (C.6). We can indeed write it as 
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so it amounts to use (C.6) with two different pole prescriptions. It is then straightforward 
to show that 



J ^((^ + fe)2 + io)(((? + fe)2-io) r(i-e)r(i-6) r(i-2e) 

2F1 (2,^-6,1- 2e, 4 V(g-")fa-")(fe- !!)(^ 



{yj{q-n){k-n)+^{q-n){k-n)y 



iVil ■ n){k ■ n) + s/^q ■ n){k ■ n)) 

(C.9) 



In the last part of the present section we report some identities for the hypergeometric 
functions appearing in the calculation: 



2 



2F1 (1, 1 - e, 1 - 2e, = (1 + zf2F^ (l, l + e, l - e, z") , \z\ < 1; 

(C.IO) 

r(l-e) , 2/^(1 - 

(C.ll) 

(^2,^-6,1- 2e, 4^^^^) = (1 + (2, 2 + e, 1 - e, z^) , \z\ < 1; 

(C.12) 

(C.13) 

To prove (C.13) we consider the following relation 

2Fi{2,2 + e,l-e,z) = (^2Fi{2,l + e,l - e, z) + z^^^2Fi{2,2 + e,2 - e, z)j , 

(C.14) 

and then we use the integral reppresentation of the hypergeometric functions in the right 
hand side of the previous equation getting 

1+. r(2-.) ,f\,y'*'(^-yr!-^\ (C.15) 



(0.16) 



l-£r(2 + £)r(-2£) J„ " (l-yzf 

Using the relation 

1 + e r(2-e) _ r(l-e) 
l-er(2 + e)r(-2e) ~ r(l + e)r(-2e) ' 

we prove (C.13). 
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